import pandas as pd
import numpy as np
import math
import matplotlib.pyplot as plt

def readDataSet(filename, div=2):
    data = np.loadtxt(filename)
    return data

def readDataSetwithLabel(filename, div=2):
    matrix = np.loadtxt(filename)
    label = matrix[:, -1]
    matrix = np.delete(matrix, -1, axis=1)
    return matrix, label

def readLabel(filename):
    data = pd.read_csv(filename, header=None)
    matrix = np.array(data)
    return matrix

def get_dist(location):
    length = location.shape[0]
    dist = np.zeros((length, length))
    ll = []
    begin = 0
    while begin < length - 1:
        end = begin + 1
        while end < length:
            d = np.linalg.norm(location[begin] - location[end])
            dist[begin][end] = d
            dist[end][begin] = d
            ll.append(d)
            end = end + 1
        begin = begin + 1
    ll = np.array(ll)
    return dist, ll

def sortLL(ll):
    return np.sort(ll)

def get_dc(sortedll, p):
    postion = int(p * len(sortedll))
    dc = sortedll[postion]

def get_rho(dist, dc):
    length = dist.shape[0]
    # 计算局部密度
    rho = np.zeros((length, 1))
    # guassian kernel
    begin = 0
    while begin < length - 1:
        end = begin + 1
        while end < length:
            rho[begin] = rho[begin] + math.exp(-(dist[begin][end] / dc) ** 2)
            rho[end] = rho[end] + math.exp(-(dist[begin][end] / dc) ** 2)
            end = end + 1
        begin = begin + 1
    return rho

def get_rho_sortind(rho):
    return np.lexsort((rho))[::-1]

def get_delta(rho_sorted_ind, dist):
    length = dist.shape[0]
    delta = np.ones((length, 1)) * np.inf
    nearest = np.zeros((length, 1), dtype=int)

    # 对于局部密度最大值取所有点的最大距离
    maxDist = -1
    maxRhoIndex = rho_sorted_ind[0]
    for i in range(length):
        if dist[maxRhoIndex][i] > maxDist:
            maxDist = dist[maxRhoIndex][i]
    delta[maxRhoIndex] = maxDist
    nearest[maxRhoIndex] = maxRhoIndex

    # 对于其他点取密度比起大的点中的距离最小值
    begin = 1
    while begin < length:
        index1 = rho_sorted_ind[begin]
        end = 0
        while end < begin:
            index2 = rho_sorted_ind[end]
            if dist[index1][index2] < delta[index1]:
                delta[index1] = dist[index1][index2]
                nearest[index1] = index2
            end = end + 1
        begin = begin + 1

    return delta, nearest

def plotRhoDelta(dir, rho, delta):
    fig = plt.figure(figsize=(10, 10))
    plt.plot(rho, delta, 'o')
    plt.xlabel('rho'), plt.ylabel('delta')
    output = dir + '/decisiongraph.png'
    plt.savefig(output)
    plt.show()

def clustering(rho_, delta_, rho, delta, nearest, sorted_rho):
    length = len(rho)
    centers = np.where(rho >= rho_ and delta >= delta_)
    cluster = np.zeros((length, 1), dtype=int)
    count = 1
    for ind in sorted_rho:
        # 从最高密度开始确定簇中心
        if ind in centers[0]:
            cluster[ind] = count
            count = count + 1
        else:
            break

    for ind in sorted_rho:
        if cluster[ind] == 0:
            cluster[ind] = cluster[nearest[ind]]

    return cluster

def get_border_sets(dist, result, dc, rho):
    length = len(result)
    num = np.max(result)
    is_border_point = np.zeros((length, 1))
    rho_b = np.zeros((length + 1 , 1), dtype=float)

    begin = 0
    while begin < length - 1:
        end = begin + 1
        while end < length:
            if dist[begin][end] < dc and result[begin] != result[end]:
                # 是border point
                is_border_point[begin] = 1
                the_cluster = result[begin]
                if rho[begin] > rho_b[the_cluster]:
                    rho_b = rho[begin]
                break
            end = end + 1
        begin = begin + 1

    return rho_b

def cut_border(rho, rho_b, result):
    length = len(rho)
    for i in range(length):
        if rho[i] > rho_b[result[i]]:
            continue
        else:
            result[i] = -1

def plotResult(location, result, filename, p):
    length = len(result)
    markers = ['.', '*', '+', 'x', '^']
    colors = ['maroon', 'red', 'peru', 'gold', 'olive', 'yellowgreen', 'lawngreen', 'springgreen']
    colors = colors + ['turquoise', 'teal', 'deepskyblue', 'dodgerblue', 'royalblue', 'navy']
    colors = colors + ['slategrey', 'orchid', 'm', 'deeppink', 'crimson']
    plt.figure(figsize=(10, 10))
    plt.title("DPC p={%.2f}".format(p))
    for i in range(0, length):
        index = int(result[i])
        if index == -1:
            plt.plot(location[i][0], location[i][1], color=(0,0,0),  marker='.')
        else:
            plt.plot(location[i][0], location[i][1], color=colors[index%19], marker=markers[index%5])
    plt.xlabel('Attribute 1'), plt.ylabel('Attribute 2')
    plt.savefig(filename)
    plt.show()

def DPC(filename, with_label, output_dir, p):

    if with_label:
        matrix = readDataSet(filename)
    else:
        matrix, label = readDataSetwithLabel(filename)

    dist, ll = get_dist(matrix)
    sorted_ll_ind = sortLL(ll)
    dc = get_dc(sorted_ll_ind, p)
    rho = get_rho(dist, dc)
    sorted_rho_ind = get_rho_sortind(rho)
    delta, nearest = get_delta(sorted_rho_ind, dist)

    decision_graph = "./" + filename[:-4] + "_decision_graph_p={.2f}.png".format(p)

    plotRhoDelta(decision_graph, rho, delta)

    rho_ = input("input rho_\n")
    delta_ = input("input delta_\n")

    result = clustering(rho_, delta_, rho, delta, nearest, sorted_rho_ind)

    rho_b = get_border_sets(dist, result, dc, rho)
    cut_border(rho, rho_b, result)

    return result, rho, delta

from getRI import rand_index
from getNMI import printNMI

if __name__ == "__main__":
    p = 0.02
    filename = './'
    output_img = 'result.png'
    output_file = 'result.txt'
    output_dir = '.'
    is_with_label = 0
    result, rho, delta = DPC(filename, is_with_label, output_dir, p)

    matrix = readDataSet(filename)
    plotResult(matrix, result, output_img)
    np.savetxt(output_file, result)

    nmi = printNMI("label.txt", output_file)
    ri = rand_index("label.txt", output_file)

    np.savetxt("rho={:.2f}_delta={:.2f}_nmi={:.2f}_ri={:.2f}".format(rho, delta, nmi, ri)+output_file, result)

